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NUMERICAL  SOLUTIONS  OF  SPHERICAL  BLAST  WAVES 
Harold  L.  Brode 

The  RAND  Corporation,  Santa  Monies,  California 

ABSTRACT 

The  strong -shock,  point -scarce  solution  and  spherical  isothermal  distri¬ 
butions  cere  used  as  initial  conditions  for  a  numerical  integration  of  the 
differential  equations  of  gas  motion  In  Lagrangean  form.  The  Von  Neumann - 
Richtsyer  artificial  viscosity  was  engloyed  to  avoid  shock  discontinuities. 

The  solutions  cere  carried  tram  two  thousand  atmospheres  to  less  than  one -tenth 
atmospheres  peak  overpressure.  Results  include  overpressure,  density,  particle 
velocity,  and  position  as  functions  of  time  and  spans.  The  dynamic  pressure, 
the  positive  and  negative  impulses  of  both  dynamic  pressure  and  static  over¬ 
pressure,  positive  and  negative  durations  of  pressure  and  velocity,  and  shock 
values  of  all  quantities  are  also  described  for  various  times  and  radial  distances. 
Analytical  approximations  to  the  numerical  results  are  provided. 

I.  INTRODOCTIOK 

The  problem  of  a  spherical  blast  in  air  has  been  solved  analytically  for 
strong  shocks  in  an  ideal  gas  by  J.  von  Neumann^,  and  In  the  weak  shock  approxi¬ 
mation  by  H.  A.  Bethe1,  For  intermediate  regions  it  has  been  found  necessary  to 

i 

resort  to  numerical  methods.  The  availability  of  high  speed  computing  machines 
has  mate  possible  the  solution  of  this  problem  in  the  Intermediate  range  with 
considerable  accuracy. 

Two  other  groups  are  currently  engaged  In  performing  this  work  by  independent 

methods.  Under  J.  Yon  Neumann  and  H.  Goldstlne  at  Princeton  a  shock  fitting 

o 

method  has  been  ceployed*  The  differencing  method  of  Peter  Lax  is  being  used 
at  Kev  York  University  by  8*  Lovell. 
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The  method  used  In  tills  paper  Is  due  to  von  Neumann  and  Richtsyer'*  and 

eoployB  an  artificial  viscosity  as  a  mechanism  for  avoiding  shock -front  Ais- 

k 

continuities.  Previously  T.  3.  Walton  reported  sons  results  using  this  method 

and  an  initial  Isobar ie  sphere  of  about  13  atmospheres. 

The  integrating  process  consists  of  the  step -vise  solution  of  difference 

agnations  which  approximate  the  differential  equations  of  motion  of  the  gas. 

In  order  for  such  a  procedure  to  be  workable,  however,  a  number  of  practical 

conditions  must  be  satisfied.  The  differencing  scheme  must  he  stable,  must 

offer  reasonably  detailed  results,  must  conserve  numerical  significance,  and  when 

put  In  the  form  of  coded  instruations  for  a  high-speed  computer,  must  be  fast 

enough  to  reach  desired  solutions  with  a  reasonable  expenditure  of  time* 

In  addition,  sealing  problems  arise  from  the  wide  range  of  numerical  values 

that  must  be  accommodated.  Furthermore,  the  continual  increase  of  shock  radius 

requires  some  recurring  adjustments  in  sene  number  and  sizes  to  preserve  a 

maximum  amount  of  Information  in  a  limited  number  of  semes.  In  the  end,  several 

eospromlses  must  be  worked  between  such  conflicting  requirements  as  the  seed  for 

frequent  printed  or  card -punched  records  of  the  progress  of  the  problem  and  a 

desire  to  minimise  the  machine  time  involved  in  printing  or  punching.  Similarly, 

the  desire  for  small  tones,  sharp  shock  fronts  and  smooth  results  is  inconsistent 

with  the  need  for  large  time  steps  and  limited  space  points.  For  simplicity, 

thesa  first  problems  were  done  with  an  homogeneous  ideal  gas  (constant  specific 

heat  ratio,  7  »  1.4),  and  for  the  simple  initial  conditions  of  either  the  strong 

1 

shock  point  source  solution  or  an  isothermal  sphere* 
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IX.  MBTHOD  FOR  HUMTOICAL  INTEGRATION 
BQUATI0K8  CBF  MOTI0H 

The  Lagrangean  equations  of  action  are  reduced  to  dimensionless  parameters 
wherein  pressure  (?) ,  density  (p) ,  and  velocity  (u)  are  measured  in  units  of 
ambient  pressure  (?.  ),  density  (p,),  and  sound  velocity  (Oj,  respectively.  In 

o  o  p 

the  following  discussion,  the  pressure  will  frequently  be  expressed  in  atmospheres 
where  one  atmosphere  is  defined  as  equal  to  the  pre -shock  ash  lent  pressure  (PQ) . 
Where  the  expression  overpressure  (AF)  is  used  the  reference  Is  to  the  pressure 
in  atmospheres  is  excess  of  tbs  ambient  pressure  (?o),  i.e.,  AP  »  P-1.  We  shall 
continue  to  apeak  of  a  velocity  (u) ,  however,  when  more  properly  we  might  refer 
to  a  aaefa  aunfaer  (u).  The  radial  dlatance  r(rc,t).  Is  expressed  in  energy- 
reduced  dimensionless  unite  (rQ  being  Lagrangean  distance,  and  t  the  tine),  such 
that 

X  »  r/«  and  *  ro/«  ,  (l) 

where  e  is  a  length  expressing  the  energy  and  ambient  pressure  sealing: 

J  ♦  f-* ,s  - IfPrir •  <*> 

°  °  o 

Etot  is  total  blast  energy  and  is  the  specific  internal  energy.  The 
subtracted  term  represents  the  pre-shoek  intemel  energy  of  the  gas  engulfed  by 
the  shook,  and  R  is  the  shock  radius.  Tina  (t)  is  defined  in  dinensionlees  units 
(*)  such  that 

t  •  t  e0/e  . 

The  artificial  viscosity  (q),  which  acts  like  a  pressure,  is  in  units  of  the 
pre-shoek  ambient  pressure  (PQ). 
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In  these  unite  the  Lagrangean  equations  of  motion  are  written  as  follows j 

"  p?  ^  It  *  "  *  |f^tp  (*“»)■  *  (3) 

I?  *  ^F~  3x  (p  +  <)  (momentum)  ,  ( k ) 


|f  *  p  It  D*  *  (y  -  1)  i]  (energy)  , 


(5) 


(6) 


In  these  equations  the  Lagrange  an  variable  Is  taken  to  be  x  -  1/3  (f^/e)^. 

&  Eq.  5  (energy  conservation)  the  internal  energy  for  an  Heal  gas  has 
been  used. 


P  po 

pXT^T 


(7) 


FORM  OP  ARTIFICIAL  VISCOSITY 

An  appropriate  viscosity  for  the  ease  of  an  outward  moving  spherical  shock 
wave  is  the  following: 


where  Ax  is  the  grid  else  and  M  is  the  number  of  grid  zones  In  the  shook  front. 

The  form  of  q  chosen  here  is  basically  similar  to  that  of  von  Neumann  and  Rlehtmyer, 
and  satisfies  similar  conditions  fen*  the  spherical  earn.  However,  there  Is  no 
convenient  steady-state  solution  by  which  to  demonstrate  the  connection  of  variables 
across  the  shook,  and  the  form  chosen  here  Is  only  asymptotically  the  earns  as  a 
verified  form  for  a  plane  problem. 
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Tbs  fact  that  the  plana  fora  approximation  la  good  for  shook  thicknesses 
which  are  small  relative  to  the  shook  radius  is  borne  out  bp  the  close  agree¬ 
ment  with  a  similarity  solution  (with  q)  investigated  by  R.  Latter^  la  which 
he  solves  numerically  the  resulting  ordinary  differential  equations  for  cases 
of  interest  to  this  problem. 

The  particular  form  chosen  for  the  viscosity  has  the  edvantege  that  it 
contributes  nothing  la  regions  of  expansion,  and  is  non -sere  only  in  the  compres¬ 
sion  phase  of  a  chock.  In  Lcgrmngean  coordinates  this  has  the  advantage  of 
eliminating  a  spurious  contribution  mar  the  origin  where  the  positive  velocity 
gradient  Is  large. 

The  artificial  viscosity,  being  essentially  a  diffusive  mechanism,  may 
force  one  to  eee  excessively  small  time  increments  in  order  to  achieve  smooth 
results  on  a  fine  space  grid. 


INITIAL  COBDmtW 

Two  general  types  of  Initial  conditions  were  taken,  (l)  a  point-source, 
and  (2)  an  Isotherm!  sphere .  A  point-source  ease  wae  run  in  the  greatest  length 
using  the  von  Neumann  strong-shock  solution  beginning  at  1600  atmos  shock  over¬ 
pressure,  and  running  down  to  less  than  0.06  atmos.  The  point-source  solution, 
starting  at  199  atmos  peak  overpressure,  was  used  to  start  the  cos  problem  carried 
out  on  the  CRDVAC  machine  at  the  Ballistics  Research  Laboratories,  Aberdeen, 
Maryland.  The  letter  problem  was  run  to  nearly  0.1  atmos  shock  overpressure. 

Tha  point-source  solution  was  also  used  ft>r  starts  at  A73  and  8l8  atmos  shook 
overpressure,  end  these  problems  were  run  down  to  lees  than  100  atmos  peek 
overpressure* 

Three  Isothermal -sphere  type  problems  have  been  run  so  far.  Two  of  these 
begin  with  a  hot  Isothermal  sphere  for  which  the  density  is  the  same  inside  and 
out.  (tee  began  at  2002  atmos,  the  other  at  121  atmos,  overpressure. 
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One  cold  Isothermal  sphere  was  also  used,  beginning  with  121  atmos  over¬ 
pressure  and  a  temperature  inside  equal  to  the  outside  temperature . 

These  three  problems  were  carried  down  to  shock  overpressures  less  than 
unity. 

Numerous  shorter  problems  were  also  run,  principally  for  the  sake  of  empirical 
checks  on  computational  methods,  or  in  amplification  of  special  features  of  the 
previously  mentioned  problems. 

Short  plane-geometry  problems  were  found  useful  in  checking  the  validity  of 

various  computational  tricks,  and  in  studying  the  effect  of  certain  changes  in 

the  difference  equations.  The  effect  of  the  size  of  the  viscosity  constant  was 

also  more  easily  studied  in  the  plane  case.  Time  increments  smaller  than  that 

i 

prescribed  by  the  Courant  Condition  (Ar  =  (pp)*^Ar)  by  a  factor  of  2  to  4  were 
found  necessary  for  smooth  results.  In  an  attempt  to  reduce  the  required  computing 
time  the  stable  differencing  scheme  of  Du  Fort  and  Frankel  for  diffusion-type 
equations  was  employed  in  place  of  explicitly  carrying  a  viscosity  quantity  (q). 
Unfortunately  such  a  scheme  has  some  practical  disadvantages.  It  requires  carry¬ 
ing  throughout  a  machine  calculation  sets  of  data  for  all  space  points  for  two 
different  times.  Furthermore,  computing,  changing  time  increments,  and  combining 
space  points  all  become  more  tedious.  Besides  these  disadvantages,  additional 
terns  must  be  introduced  to  correct  for  the  excess  energy  introduced  by  the 
differencing  scheme.  On  the  other  hand,  the  very  general  nature  of  the  viscosity 
method,  the  ease  of  its  applicability,  and  the  precision  with  which  it  reproduces 
the  Hugoniot  conditions  across  a  shock  would  seem  to  offset  the  more  stringent 
time  requirement.  Use  of  this  method  for  non-ideal  gases  is  not  considered  in 
this  paper,  however. 

It  is  frequently  convenient  to  use  unequal  zone  sizes.  For  instance,  the 
use  of  small  zones  through  the  shock  front  provides  a  sharp  shock  at  very  little 
cost  in  computing  time.  The  use  of  such  unequal  zones  was  empirically  validated 
in  this  problem  by  repeating  calculations  with  quite  different  zone  choices. 
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The  else  of  the  tin  increments  were  automatically  doubled  whenever  the 
stability  conditions  would  allow  It.  Two  conditions  exist,  one  being  the  usual 
Courant  Condition  and  the  other  a  diffusion-type  condition  loosed  by  the  artificial 
viscosity 


At  ^  (Ax)2  -|- 


(9) 


The  total  energy  in  the  blast  wave  must  be  conserved,  but  a  check  on  the 
total  energy  is  not  a  very  sensitive  test  of  the  correctness  of  the  results. 

Several  machine  errors  were  detected  (and  eliminated)  which  did  not  substantially 
effect  the  total  energy.  The  most  reliable  check  seems  to  be  to  re-run  a  problem 
with  minor  alterations.  Such  a  re-run  not  only  catches  machine  errors  but  serves 
to  measure  the  Ices  of  significance.  Unfortunately  this  procedure  is  costly  of 
machine  time,  end  a  Compromise  might  he  arranged  whereby  the  machine  is  instructed 
to  double  compute  and  check  itself  every  cycle.  Double  computing  cannot  eliminate 
errors  in  storing  or  transferring,  however,  and  on  the  particular  problem  reported 
here  the  system  employed  was  to  back  off  and  re-run  wherever  suspicious  aberrations 
occurred  and  to  re-run  and  ea&pare  the  entire  problem  with  a  different  zone  spacing, 
different  Viscosity,  different  time  increments,  as  well  as  scorn  differences  in 
the  difference  equations. 


The  differential  equations  are  approximated  by  the  following  difference 
equations! 

>2 


-  *  *  ♦ 


-  rU +  V 


(10) 
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III*  RESULTS 


PEAK  0VEHPRE88UBE  (Pj  VERSUS  SHOCK  RADIUS  (X  ) 

.8  ...'',  _ _ _ 8 


The  strong-shoek  point -source  solution  provides  that  the  shock  overpressure 
(AP8)  should  depend  oa  the  inverse  eube  of  the  shook  radius  (X^).  The  problems 
that  sere  begun  with  strong-shock,  point -source  values  continued  to  obey  the 
inverse  cube  law  down  to  10  atmes  at  which  point  the  overpressure  had  become  3 
percent  higher  than  the  strong  shock  prediction.  The  addition  of  1  atmos  to  an 
inverse  cube  term  gives  a  relation  valid  over  a  greater  range  of  overpressures, 
deviating  from  the  calculated  eurve  by  lees  than  5  percent  at  5  atmos. 
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AP  *  O.I567  +  1  atmos  (16) 

0  8 

The  shock  radius  fc  )  is  la  the  dimensionless  units  ( energy /pressure -reduced ) 
described  ia  Sect lea  II. 

At  still  loner  pressures  the  following  enpirieal  fit  (applies. 

Af  -  Sa|21  +  2aM2  +  -  0.019  atmos  ,  (17) 

K  K  8 

for  0.1  <  AP  <  10,  or  0.26  <  V  <  2.8  . 

For  all  the  calculations  begun  with  a  strong-shock  the  solid  curve  in 
Fig.  1  presents  the  peak  overpressure  as  a  function  of  shock  radius. 

The  dashed  curves  represent  the  peak  overpressures  resulting  frost  the  hot 
isothermal  spheres  with  normal  density  inside  (beginning  at  2002  and  122  atmos 
overpressure) .  Hot  ice  that  the  overpressure  becomes  Indistinguishable  from 
that  of  the  strong  shock  at  a  radius  where  the  mass  of  air  engulfed  by  the  shock 
Is  10  times  the  initial  mass.  In  the  strong  shock  region  this  corresponds  to  a 
drop  in  overpressure  to  about  17  percent  of  the  initial  value,  and  in  any  case 
to  a  radius  a  little  more  than  twice  the  initial  radius. 

The  dotted  curve  represents  the  overpressure  from  an  initially  cold  iso¬ 
thermal  sphere  (nonail  temperature,  high  pressure) .  Because  of  the  slower 
rarefaction  speed  the  shock  pressure  does  not  rise  to  the  point-source  value 
until  quite  late,  l.e* ,  until  the  rarefaction  which  moves  inward  frees  the 
initial  pressure  front  has  reached  the  center. 

In  each  case* of  an  initial  isothermal  sphere  there  are  some  slight  oscil¬ 
lations:  About  the  point  source  overpressure  curve  which  arc  caused  both  by 
rarefactions  and  small  shocks  which  foam  behind  the  front  shock,  move  inward 
and  reflect  off  the  origin  and  then  move  out  to  overtake  the  shock  front. 

(See  a  later  discussion  of  these  special  features  of  the  sphere  problem.) 
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PEAK  DYNAMIC  PRESSURE  (Q^)  AND  PEAK  PARTICLE  VELOCITY  («  )  VERSUS  SHOCK  RADIUS  (A  ) 


This  single  form  is  reflected  in  tbs  peak  dynamic  pressure  (Q  »  l/a  p  u^) 
which  falls  somewhat  faster  than  an  inverse  cube  law  deem  to  equally  lew  shock 
pressures  (Pig*  3)  •  The  sloe  decrease  in  shook  density  (frees  6  for  strong  shocks 
to  one  for  weak  shocks)  Is  responsible  for  the  decay  being  steeper  than  the  inverse 
cube  of  the  shock  radius. 


PRESSURE  (P)  VERSUS  UGRABGEAN  DISTANCE  (Rq) 

The  pressure  behind  the  shock  wave,  shown  in  Fig.  k  as  a  function  of  the 
lagrangtan  (or  Initial  position)  variable,  retains  the  strong-shock  form  until 
quite  low  pressures.  Note,  for  exaaple,  that  the  ratio  of  central  pressure  to 
shock  pressure  remains  37  percent  down  to  20  atmoe  and  decreases  slowly  below 
that  to  33  percent  by  3  atmoe.  Beyond  3  atmoe  a  negative  phase  develops  with 
the  pressure  falling  to  as  low  as  0.8  atmoe  near  the  center. 

In  Fig.  h  the  lagrangean  position  is  given  in  arbitrary  units  (Rfi), 
related  to  the  dimensionless  unit  (X0)  by  a  constant  multiplier. 


«x 


,  V 
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(20) 


It  may  help  In  visualizing  dime ns Ions  to  assume  R^  to  be  in  centimeters;  then 

the  correspoading  blast  energy  rill  be  about  that  of  200  lbs.  of  TOT  (4.2  x  lO^erga), 


PARTICLE  VELOCITY  (tt)  AUD  DKHSITY  (p)  VERSUS  LAGRAJTGEAN  DISTAHCE  (Rq) 

Figs.  5  and  6  indicate  the  progression  of  particle  velocity  and  density 
as  functions  of  the  Lagrange an  position.  Again  the  strong-shock  fora  is  dominant 
until  as  lor  as  3  atnos  of  steak  overpressure. 

The  particle  velocity  transforms  gradually  from  its  very  linear  fora  to  one 
much  like  the  overpressure  at  large  distances  as  the  shock  rave  goes  from  strong 
to  weak. 

The  density,  which  in  the  strong  shock  is  zero  at  the  origin,  implying  an 
infinite  temperature ,  remains  zero  there  sines  no  msehaniem  ifi  included  (no 
conduction  or  radiation)  far  dissipating  this  high  temperature. 

The  density  profiles  in  Fig.  6  represent  strong  shook  initial  conditions 
at  200  atoos  steak  pressure*  The  dip  in  the  curves  that  sit  at  the  same  mess 
point  (R0  -  130)  is  due  to  the  sudden  inclusion  of  a  finite  atmosphere  ahead 
of  the  shock,  l.e.,  due  to  the  abandonment  of  the  etrong-shoek  (infinite  pressure 
ratio)  aasmaption  at  that  point.  The  Hugoniot  relatione  give  the  chock  density 
In  terms  of  the  shock  pressure  as 

p  •  rrlfor  7  m  x*k  *  *al> 

For  a  strong  shock  (?  «  so)  the  density  ratio  is  six  (p  •  6),  but  at  two 
hundred  atmoe  (P  *  200)  the  ratio  is  only  5.83. 
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71m  strong  shock  solution  is  derive!  on  the  assuBptlo&  that  the  atmosphere 
ahead  of  the  shock  wave  has  negligible  effect*  The  effect,  when  non-negligible , 

Is  to  raise  the  tesperature  through  the  shook  to  a  higher  value  than  that  given 
by  the  strong  shook  solution,  i.e*,  e  finite  shock  is  hotter  than  would  be 
predicted  by  the  "strong -shock"  theory. 

GAS  VAKEABLB8  AS  A  FUBCTIOH  OF  HCTlERXAN  POflXTIOH  (l) 

The  pressure,  density,  compression,  and  particle  velocity  ere  shown  relative 
to  their  peak  or  shock  values  as  functions  of  the  Eulerlan  position  (1)  In  Fig.  7. 
The  strong  shock  form  dominates  the  first  two  sets,  while  the  later  ones  show 
the  characteristic  positive  phase  followed  by  a  longer,  weaker  negative  phase 
and  eventually  by  a  return  to  near  pre -shock  values  at  the  origin. 


POSITIVE  ABB  SEOATIVE  PHASES  VERSUS  DI3TAHCB 

The  durations  of  the  positive  phases  for  pressure  and  particle  velocity  are 
shown  in  Fig.  6.  Although  these  durations  should  approach  the  same  value  at 
large  distances,  they  still  differ  by  7  percent  at  a  distance  of  A  »  3.0  (£?  *  0.09). 

o 

The  figure  Indicates  only  the  values  for  the  point  source  solution. 

The  duration  of  the  negative  pressure  phase  (D~) ,  unlike  that  of  the  positive 

Mr 

phase,  is  nearly  Independent  of  distance,  and  has  an  average  value  of  1.22.  The 
negative  durations  suffer  somewhat  from  loss  of  nuncrleal  significance  at  late 
times  in  the  calculation. 


TIMS  BEPEKBKHCB  OF  PH^SOSS 

In  Fig.  9,  the  curves  of  pressure  versus  time  at  various  distances  ere 
given  in  units  of  the  peak  overpressure  and  the  positive  duration  In  order  to 
Illustrate  the  change  in  the  rate  of  decay  behind  the  shock.  At  Increasing 
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distances  the  drop  la  pressure  Is  both  relatively  and  absolutely  slower  since 
the  positive  phase  duration  is  also  Increasing.  Again,  these  curves  are  for 
the  point -source  solution. 

An  exponential  form  approximates  the  time  dependence  of  the  pressure  pulses 

H*  »  (1  -  *)«-*  {22) 

8 

Here  %  is  the  time  after  shook  arrival  in  units  of  positive  duration,  APg  is 
the  shook  overpressure,  and  a  is  independent  of  Z*  This  form  is  satisfactory 
for  overpressures  less  than  one  atmosphere  {X  >  0.7*0  >  with  the  eoefflelent  (a) 
specified  by 

a  #  1/2  +  AP  ,  AP  ^  1  (23) 

For  shock  overpressures  greater  than  one  atmosphere  the  decay  is  not  a 
simple  exponential,  since  the  early  portion  requires  a  larger  a  than  the  later 
part.  Allowing  the  coefficient  a  to  be  a  function  of  the  time  (z) ,  it  may  be 
approximated  by 

a  -  1/2  +  APg  [l.l  *  (*13+  *2  APg)z]  (2k) 

for  AP  less  than  3  atmos,  and  for  overpressures  from  3  up  to  50  atmos  by  the 

8 

form 

°  *  a  +  &&  (25) 

where 

a  »  f-  0.231  +  0.388  APg-  0.0332  AP2  .for  AT  4  10 
(  0  for  APg  >10 

b  *  Ta|(0.88  +  0.(372  A|)  for  Ap^  <  10 
[a|(1.67  -  0.011  A|)  .far  Alt  ^  10 
e  -  8.71  +  0.18X3  -  10X/(APg  +  10) 


RM-1363-AEC 

9-29-54 

*  14  - 

The  tine  dope nde see  of  the  pressure  la  the  negative  phase  may  be  approximated 
by  the  fora 

AP  «  l4  AP  <0(1  -  »)€ “^  for  0.1  <  AP  <  200  (26) 

*  B 

where  AP.  Is  the  peak  negative  overpressure  (see  Fig.  10),  and  where  <*>  Is  the 
time  measured  from  the  end  of  the  positive  phase  la  units  of  the  negative  phase 
duration  (B‘).  *  ®or*  •«“•*  m  uould  .now  th.  «*»»t  to  *«*.«.  to  «ro 
as  the  shock  strength  goes  to  zero. 


DYHAMIC  PRESSURE  VERSUS  TIMS 

o 

Of  Interest  also  Is  the  fora  of  the  dynamic  pressure  (Q  *  1/2  pu  ) .  Time 
plots  of  this  function  appear  la  Fig.  11.  R»se  eurves  are  also  normalised, 
but  with  the  positive  duration  of  the  particle  velocity  (D*)  (Fig.  8),  and  the 
peak  value  of  the  dynamic  pressure  (Q  )  (Fig.  3)  • 

B 

A  similar  exponential  fora  approximates  the  dynamic  pressure, 

§--  U  -  «■*  ■  <S7> 

where  p  is  indspeadent  of  the  time  (8)  for  shocks  of  lese  than  one  atmos  peak 
overpressure  (AP.  <  l).  The  coefficient  0  changes  with  the  shock  strength, 
however: 

^  -  0.75  +  3*2  AP.  AP.  «  1  .  (28) 

B  8 


Where  the  shock  is  stronger,  a  modification  similar  to  that  given  in 
Eq.  (23)  for  the  overpressure  is  appropriate: 

p  »  d  +  f 


(29) 


for 


d  • 


50  >  AP.  >  1  , 

s 

(•  1-33  APS  fer  APfl«  3 
•  5*6  +  0.63  APa  £r  5<  AP8  <:  10 


for 


AP.  >  10 

B 


and 
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f  =  6.1|0  AP 

S 

g  -  0.7#  &PS 

These  approximate  forms  agree  with  the  numerical  values  of  Q  and  &P  in 
the  positive  phase  to  within  10  percent  (and  over  most  of  the  range  to  less 
than  L  percent)  for  values  of  APg  less  than  ten  atmos. 


POSITIVE  IMPULSE 


The  integrated  positive  overpressure  (I  )  and  the  total  positive  drag 

P 

pressure  (I*)  decrease  with  distance  from  the  blast  source  in  the  manner  shown 


in  Fig.  1£* 


u 

<;-r 


AP(t)  dt  , 


(30) 


♦  _  1 


u 


pu^  dt 


(Because  of  the  units  employed  in  this  paper*  Iu  dpes  not  become  the  usual 


dynamic  impulse  until  it  ia  multiplied  by 


Pn* 


oO 


) 


Both  of  these  impulses  may  be  fitted  by  simple  powers  of  the  radial 
distance  for  shock  overpressures  less  than  two  atmospheres? 


I+  =  0.01*3  X_1  ,  AP  <  £ 
P  s 


(3D 


i+  =  0.00U  x"2,7 

u 


NEGATIVE  IMPULSE 

The  negative  overpressure  impulse  in  the  range  below  £0  atmos  of  peak 
overpressure  can  be  expressed  to  within  percent  by 


I 

P 


r  &P  D" 
c  -  P 


(3£) 
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Since  the  duration  of  the  negative  pressure  is  constant  at  nj'  =  1.22, 

1^  »  0.61  AP_  .  (33) 

For  distances  greater  that  X  =  0.5  the  negative  peak  overpressure 
approaches  zero  inversely  as  the  distance 

AP_  =  -  0.086  X"1  ,  X  >  0.5  ,  (34) 

and  the  negative  pressure  impulse  goes  approximately  like 

i;  =  0.052  X"1  ,  X  >  0.5  .  (35) 

SPECIAL  FEATURES  OF  THE  ISOTHERMAL  SPHERE  PROBLEMS 

In  the  gas  dynamics  resulting  from  the  release  of  an  initially  static 

high-pressure  sphere,  an  inward  moving  shock  forms  behind  the  rarefaction 

wave  that  first  runs  in  from  the  surface  of  the  sphere .  This  inward- 

directed  shock  was  predicted  by  Wecken^  and  discussed  by  MeFadden®,  Shardin^ 
k 

add  others  .  It  does  not  acquire  a  net  inward  velocity  until  the  rarefac¬ 
tion  has  reached  the  center,  but  after  that  it  moves  in  and  reflects  at  the 
origin,  and  then  races  outward  to  eventually  overtake  the  main  shock. 

This  second  shock  grows  from  zero  strength  to  presumably  an  infinite 
pressure  ratio  at  the  origin.  On  reflection,  this  shock  moves  outward, 
decaying  in  strength  about  as  the  Inverse  first  power  of  its  distance  from 
the  origin.  The  general  expansion  and  outward  motion  of  the  gas  behind 
the  main  shock  is  responsible  for  this  apparently  modest  decay  rate. 

When  this  second  shock  encounters  the  outer  surface  of  the  gas  that 
was  initially  inside  the  sphere,  a  transmitted  shock  continues  out  and  a  re¬ 
flected  shock  is  sent  inward.  The  transmitted  shock  overtakes  the  main 
shock  and  increases  it  by  as  much  as  20$.  Progressively  weaker  shocks 
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follow  from  the  reflections  at  the  Interface  and  again  at  the  origin.  This 
repeated  shocking  of  the  gas  near  the  origin  changes  the  temperature  profile 
at  later  times  from  its  initial  Isothermal  nature  to  nearer  the  point-source 
distribution  which  is  characterized  by  a  high  central  temperature  falling 
off  rapidly  with  radius. 

Some  detail  is  lost  in  the  vicinity  of  the  origin  since  the  shocks  are 
always  spread  over  a  number  of  mesh  points.  In  one  case  the  problem  was  re¬ 
run  with  mesh  sizes  about  one  fourth  the  original  size,  and  the  inward-moving 
shock  showed  some  appreciable  discrepancy  near  the  origin.  However,  this 
discrepancy  may  be  attributed  mainly  to  the  difficulty  in  identifying  the 
shock  when  the  rounding  is  comparable  to  the  shock  radius.  After  reflection 
the  difference  vanishes  again  in  spite  of  the  quite  different  histories  near 
the  center. 

IV.  CONCLUSIONS 

POIHT  SOURCE  SOLUTION 

At  increasing  distances  from  a  finite  but  sudden  source  of  energy,  the 
resulting  blast  wave  will  appear  more  and  more  like  that  from  a  point  source. 

The  blast  resulting  from  an  initial  isothermal  sphere  of  gas  at  rest  will 
assume  the  general  shape  and  values  of  the  point  source  solution  (to  within 
10  percent)  after  the  shoek  wave  has  engulfed  a  mass  of  air  10  times  the 
initial  mass  of  the  sphere.  Prior  to  this,  the  shock  strength  is  less  than 
that  of  the  point  source  shock,  and  the  inward-traveling  rarefaction  has  not 


reached  the  center 
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A  point  source  should  leave  a  higher  temperature  and  consequently  a  larger 
percentage  of  energy  near  the  origin,  fills  energy,  no  longer  available  to  the 
shock  wave,  should  effect  a  reduction  in  the  shock  radius  (for  a  given  over* 
pressure).  However,  no  appreciable  difference  in  the  low-end  of  the  shock 
overpressure -radius  relation  appears  on  comparing  the  point  source  and  isothermal 
sphere  solutions.  In  fact,  perhaps  due  to  multiple  stocking  of  the  inner  regions 
in  the  isothermal  sphere  problems,  the  distribution  of  residual  energies  (per 
unit  volume)  and  pressures  axe  nearly  identical  around  the  origin  at  a  time 
when  the  shock  has  progressed  to  6  or  7  times  the  initial  radius*  Although 
temperatures  will  remain  different,  since  the  point  source  has  an  infinite 
temperature  at  the  origin  (and  zero  density),  the  average  temperature  (or  density) 
of  that  gas  initially  inside  the  isothermal  sphere  approaches  (within  10  percent) 
the  average  temperature  for  a  corresponding  mass  around  the  point  source. 

For  a  source  of  considerable  initial  mass,  the  peak  pressure  may  become 
quite  small  before  the  shock  has  angulfed  a  mass  of  gas  10  times  larger,  and 
the  blast  wave  may,  therefore,  remain  quite  different  from  the  point  source 
solution  throughout  regions  of  interest,  flush  is  the  case  to  scam  extent  with 
high  explosives  where  Initial  charge  shapes  will  influence  the  blast  wave  at 
all  significant  pressures.  It  is  true  to  an  even  greater  extent  in  the  spherical 
equivalent  of  a  shock -tube  type  of  blast,  where  a  gas  at  high  pressure  but  at 
normal  tesperature  is  suddenly  released,  as  in  the  normal  temperature  isobsrle 
sphere  problem  (dotted  curve  Fig.  1)  described  in  flection  I. 

APPLICATION  TO  BLAST  IN  AIR 

The  Ideal  gas  assumption  is  reasonably  valid  in  air  for  Chock  pressures 
less  than  10  atmes.  Above  that  the  gamma  in  the  expression  feu*  the  internal 
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energy  (Bq.  7}  ranges  down  to  as  low  as  1.133  end  up  to  the  mono-atomic  value 

of  1.667. 

But  slnea  the  mass  of  gae  which  has  experienced  shocks  stronger  than  10 
atffios  Is  a  small  part  (5  percent)  of  the  engulfed  mass  by  the  time  the  shook 
overpressure  is  down  to  1  atmoe,  the  blast  wave  for  constant  gamma  should  be 
reasonably  correct  for  air  at  most  Interesting  shock  pressures. 

In  the  ease  of  shook  wave  problems  Involving  only  one  space  variable, 
the  artificial  viscosity  technique  for  numerical  Integration  appears  to  be  very 
satisfactory  over  large  ranges  of  pressure  and  entropy  change.  Unfortunate ly , 
it  cannot  be  expected  to  yield  details  of  a  shock  wave  l wringing  on  a  singular 
point  such  as  the  origin  la  spherical  or  cylindrical  geometry,  since  the  nature 
of  this  method  Is  such  as  to  spread  the  shook  front  over  a  number  of  grid  points* 
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FIGURE  CAPTIOUS 

-  Radial  dependence  of  the  peak  of  shock  overpressure.  The  solid 
curve  represents  the  point  source  solution.  The  dashed  curves 
represent  the  restate  from  two  different  initial  isothermal  spheres, 
one  at  8000  atmos  and  the  other  at  121  atmos,  both  spheres  were  hot 
with  equal  density  inside  and  out.  The  dotted  curve  represents 

the  resulting  shock  overpressure  from  an  initially  cold  isobar ic 

sphere  with  tenperature  equal  inside  and  out.  Tbs  pressures  are 

in  atmos  and  the  distances  in  units  of  (E.  . /p 

tot  o 

•  Particle  velocity  or  Mach  number  at  the  shock  as  a  function  of 
shock  radius  for  the  point  source  solution.  Velocity  is  in  units 
of  the  pre -shock  sonic  velocity  and  tha  radius  is  In  unit®  ,  of 

<WV1/3- 

p 

-  Peak  or  shock  dynamic  pressure  (Q  *  1/2  pu  )  versus  shock  radius  for 

the  point  source  solution.  The  dynamic  pressure  is  in  atmos  and  the 

radius  in  units  of  (E.  . /pj1//3. 

tot  o 

-  Pressure  in  atmos  as  a  function  of  the  Lagrange  or  mass  position 
for  the  point -source  solution  at  times  indicated.  The  position  is 
in  arbitrary  units  of  (Etot/Po)1^/l627.2,  and  the  time  la  in  units 
©f  <WPo>1/3/Co  ll!Mm  Co  18  the  P*«*hoek  sonic  velocity, 

-  Particle  velocity  in  units  of  preshock  sonic  velocity  as  a  function 

of  Lagrange  or  mass  position  for  tbs  point-source  solution  at  the 

times  indicated.  The  position  Is  In  units 

and  the  time  is  in  units  of  (E.  ./¥  . 

tot  o  © 
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Fig*  6  -  Density  in  units  of  preshock  density  as  a  function  of  Lagrange  or 
mas  position  for  the  point-source  solution  at  the  tints  indicated. 
The  position  is  in  units  of  (Etcrfc/P0) 


Fig.  t  *  Hydrodynamic  quantities  in  units  of  their  peak  values  as  a  function 
of  Sulerlan  or  space  position  at  the  times  indicated.  The  radius  is 
in  units  of  The  solid  curves  represent  overpressure 

(&p/aP  ),  the  long  dash  curves  represent  particle  velocity  (u/u  ), 
the  dot-dash  curves  represent  density  (p/p„),  and  the  short  dash 
curves  represent  coapression  [(p-l)/(p8*l)| *  The  shook  values  of 
the  overpressure  (£P  )  for  these  tines  are  121.5*  20.10,  2.03,  1.01, 
0.33&*  0.0701  ordered  on  increasing  tins.  For  the  particle  velocity 
the  peak  values  are  0.k9,  3*37*  O.873,  0.52k,  0.216,  0.0525,  and  for 
the  density  the  peak  values  safe  5*66,  k.68,  2.135*  1.625*  1.228, 
1.0508. 


Fig.  8  -  Duration  of  positive  phase  for  pressure  (Dp  and  particle  velocity 

<D+)  versus  distance  { Euler  ian)  where  tine  is  in  units  of  (E,  . /P)1^ 
u  tox  o 

and  distance  la  in  units  of  (Bto%/Pe)^. 

Fig.  9  *  Pressure  as  a  function  of  tine  where  the  pressure  is  in  units  of 
the  peak  pressure  (£P)(Fig,  1)  and  the  tins  (Z  )  is  in  units  of 

0  P 

the  positive  duration  (Dp  (Fig.  8).  The  numbers  indicate  the 

Jr 

corresponding  peak  overpressures  (dPg). 

Fig.  10  -  Peak  negative  overpressure  in  atnee  as  a  function  of  radial  distance 


in  units  of  (E^ot/Po) 


1/3 
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Pig,  11  -  Dynamic  pressure  (Q  «  l/2  pu2)  as  a  function  off  time  where  tins 
pressure  is  in  units  of  the  peak  dynamic  pressure  (Q^  (Fig.  3) , 
and  the  time  is  In  units  of  tint  positive  duration  of  the  velocity 
'  (Dy)  (Fig.  8).  The  humbersindieate  thecorresponding  peak  over- 
pressures  (dF^.  This  rapid  decrease  in  effective  pulse  duration 
with  increasingpeak  overpressure  has  the  effect  of  making  the 
higher  dynamic  impulses  much  t£a,  sharpest.  ' 

?ig.  12  -  Tbs  radial  dependence  of  the  positive  pressure  impulse  (ij)  the 
negative  pressure  impulse  (ij),  and  the  positive  dynamic  impulse 
Iq  in  units  of  atmospheres  (where  time  is  in  the  dimensionless 
units  (®t0^/?0)^^AS0)  the  radial  distance  is  in  units  of  (s^/p^)1^. 
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